##########################################################################################

library(data.table)
library(optparse)
library(ggplot2)
library(Seurat)
library(ggsci)
library(ArchR)
library(ggplotify)

##########################################################################################
option_list <- list(
    make_option(c("--cluster"), type = "character"),
    make_option(c("--rna_file"), type = "character"),
    make_option(c("--atac_file"), type = "character"),
    make_option(c("--cell_cluster_file"), type = "character"),
    make_option(c("--out_path"), type = "character") 
)

if(1!=1){
    
    cluster <- "cluster2"
    ## 单细胞表达文件
    rna_file <- "~/20231121_singleMuti/results/qc_atac_v3/germ/testis_combined.annotationCellType.qc.Rdata"

    ## atac文件
    atac_file <- "~/20231121_singleMuti/results/qc_atac_v3/germ/testis_combined_peak.combineRNA.qc.Rdata"

    ## sperm重新聚类的文件
    cell_cluster_file <- "/public/home/xxf2019/20231121_singleMuti/results/celltype_plot/recluster_allcell/cluster2_atac_cluster.tsv"

    ## 输出
    out_path <- "~/20231121_singleMuti/results/celltype_plot/recluster_allcell"

}

###########################################################################################
parseobj <- OptionParser(option_list=option_list, usage = "usage: Rscript %prog [options]")
opt <- parse_args(parseobj)
print(opt)

cluster <- opt$cluster
rna_file <- opt$rna_file
atac_file <- opt$atac_file
out_path <- opt$out_path
cell_cluster_file <- opt$cell_cluster_file

###########################################################################################

a <- load(rna_file)
DefaultAssay(scrnat) <- "MAGIC_RNA"
## scrnat

b <- load(atac_file)
## testis_combined_peak_combineRNA
projHeme5 <- testis_combined_peak_combineRNA

sperm_cell_cluster <- data.frame(fread(cell_cluster_file))

###########################################################################################

sperm_cell_cluster$cell_rna <- gsub( "#" , "_" , sperm_cell_cluster$cell )
#c1_cell <- subset(sperm_cell_cluster , cluster == "C1")
#c2_cell <- subset(sperm_cell_cluster , cluster == "C2")

#col_sample <- c(
#    rgb(red=212,green=31,blue=37,alpha=255,max=255) ,
#    rgb(red=39,green=45,blue=106,alpha=255,max=255) 
#    )
#names(col_sample) <- c("C1" , "C2")

###########################################################################################
## RNA
scrnat$sperm_atac_recluster <- "other"
projHeme5@cellColData$sperm_atac_recluster <- "other"

for( clus in unique(sperm_cell_cluster$cluster) ){

    clus_cell <- subset( sperm_cell_cluster , cluster == clus )
    scrnat$sperm_atac_recluster <- ifelse( scrnat$cell %in% clus_cell$cell_rna , clus , scrnat$sperm_atac_recluster )
    projHeme5@cellColData$sperm_atac_recluster <- ifelse( rownames(projHeme5@cellColData) %in% clus_cell$cell , clus , projHeme5@cellColData$sperm_atac_recluster )

}


p1 <- DimPlot(scrnat, reduction = "umap", label = TRUE, group.by = 'sperm_atac_recluster') + theme(legend.position = 'none') 
p2 <- plotEmbedding(ArchRProj = projHeme5, colorBy = "cellColData", name = "sperm_atac_recluster", embedding = "UMAP", plotAs="points" )

###########################################################################################
## 对RNA的细胞重新聚类
scrnat_sperm <- subset(scrnat , seurat_clusters==gsub( "cluster" , "" , cluster ))
scrnat_sperm <- FindNeighbors(scrnat_sperm, reduction = "pca", dims = 1:30)
scrnat_sperm <- FindClusters(scrnat_sperm , graph.name = grep( "snn" , names(scrnat_sperm@graphs) , value = T ) , resolution = 0.6)

## 重聚类
scrnat_sperm <- RunUMAP(scrnat_sperm, dims = 1:10)

p3 <- DimPlot(scrnat_sperm, reduction = "umap", label = TRUE, group.by = 'sperm_atac_recluster') + theme(legend.position = 'none') 

out_file <- paste0( out_path , "/" , cluster , ".AllCell.atac_recluster.pdf" ) 
pdf(out_file)
print(p1)
print(p2)
print(p3)
dev.off()
